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Abstract. We investigate methods to estimate the parameters of the 
gravitational-wave signal from a spinning neutron star using Fourier transformed 
segments of the strain response from an interferometric detector. Estimating the 
parameters from the power, we find generalizations of the PowerFlux method. 
Using simulated elliptically polarized signals injected into Gaussian noise, we 
apply the generalized methods to estimate the squared amplitudes of the plus 
and cross polarizations (and, in the most general case, the polarization angle), 
and test the relative detection efficiencies of the various methods. 



PACS numbers: 04.80.Nn, 95.75.-z, 07.05. Kf 

1. Introduction: parameter estimation using power 

A periodic gravitational wave incident on an interferometric detector will produce a 
strain response of the form 

h{t) = A+F+ (V>, t) cos $(t) +A X F X t) sin $(f), (1) 

where h(t) is the strain, A + and A x are the amplitudes of the plus and cross 
polarizations of the gravitational wave, F + and F x are the respective sky-position 
dependent response functions (or antenna patterns) of the detector, tf> is the 
polarization angle, and $ is the gravitational wave phase, which contains modulations 
from Doppler shifts due to the relative motion between the source and the detector 
and the frequency evolution of the source [1] . (The response of a bar detector can be 
written in a similar form, but with antenna patterns differing from those in [1].) 

The LIGO Scientific Collaboration (LSC) has used fully-coherent and semi- 
coherent methods to search for unknown sources of these waves, while also developing 
hierarchical schemes that use a combination of these methods (see [H [3] and references 
therein). The starting procedure for these methods is to divide the strain data into 
short segments and take the discrete Fourier transform of each, to create Short Fourier 
Transforms (SFTs) of the data. The time-baseline used to generate the SFTs is 
designated Tsft- This is typically chosen to be 30 min, so that signals with frequencies 
less than approximately 1000 Hz from an isolated source will not, within this time 
period, shift in frequency by more than half the width of an SFT frequency bin. 
(The method can also be applied to sources in binary systems, by using a shorter 
time-baseline dependent on the size of the expected orbital Doppler shifts.) 



Generalized PowerFlux methods ... periodic gravitational waves 



2 



The PowerFlux method 0] is a semi-coherent method for detecting periodic 
gravitational waves using a weighted average of the power from SFTs. Compared 
with the StackSlide and weighted Hough methods, that also combine power from 
SFTs, in general the PowerFlux method has the highest sensitivity |5>]- Hierarchical 
and fully-coherent searches that use a coherent time-baseline longer that 30 min are 
more sensitive, but also more computationally costly. Thus, PowerFlux remains a 
vital search method for periodic gravitational waves. 

PowerFlux is able to estimate the squared amplitude (either A\ or A 2 X ) of 
a linearly polarized signal (A x = or A + = 0) or a circularly polarized signal 
(A 2 , = A 2 X ), which can also be used as a detection statistic. However, it does not 
directly estimate the parameters of elliptically polarized signals (A\ ^ A x / 
in general), nor does it directly recover the polarization angle ip. In this paper we 
investigate whether it is possible to extract this extra information from a potential 
signal present in the SFTs and, if so, whether we can hnd a method using power from 
SFTs that is more efficient at detecting signals than standard PowerFlux. 

We have previously (in [5]) investigated using the real and imaginary parts of 
the SFTs to estimate A+, A x , and -0 using the method in pQ; the method worked 
well for searches with long coherent time-baselines, but failed to be robust on time- 
baselines of 30 min, for a number of reasons, as detailed in [5j. Here, we show instead 
how to estimate the parameters of elliptically polarized signals using the power in 
SFTs, thus generalizing the PowerFlux method. As a first step, we present a new 
derivation of the PowerFlux method in the next section. The subsequent sections 
give two generalizations, and then a comparision of the detections efficiencies of the 
methods (and the distributions of the estimated parameters). The last section gives 
another way to form a detection statistic using SFT power, analogous to finding a 
maximum likelihood statistic, and suggests future work. 



2. Derivation of the PowerFlux method 



By definition of Tsft, we can treat F+, F Xl and the frequency of the signal at the 
detector as constant over the duration of an SFT. Ignoring losses due to the difference 
between this frequency and the closest SFT bin frequency, the normalized power in 
the signal is 

|^=0.5(^+A 2 X J^)T SFT> (2) 

J SFT 

where h is the discrete Fourier transform of h(t) divided by the sample rate of the 
data, and it is understood that F + and F x are evaluated at the midpoint of each SFT. 
Equation ([2|) represents the expected signal power for an elliptically polarized signal 
from one SFT. We label the SFTs using index a, and consider searching for a linearly 
polarized signal (with A x — 0) which is present in the sequence of SFTs. Given the 
SFT data, x a , we define 

Pa = ^ (3) 
J SFT 

to be the power taken from the SFT bin closest to the expected signal frequency, and 
_ ^ [P a - 0.5AlFl a T SFT ]i 
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to be the noise-weighted sum of the square deviations in power between the signal and 
the data; S a are the one-sided power spectral densities of the noise for the frequency 
bins used in each corresponding SFT. 

A natural way to estimate A+, analogous to \ 2 minimization, is to find the value 
of A\ that minimizes g. We therefore solve 

9g_ = _J- (P ° Z 0^A 2 + F 2 a T SFT )F 2 a T SFT = q (5) 



9A\ ^ S% 



for A\, and obtain 



^ = 4E^i^/E^- (6) 

a °a J SFT a °a 

Equation ^ is the detection statistic for the PowerFlux method given in [3J [3J , 
although the derivation given here is different. We refer to this as the linear PowerFlux 
method. To instead apply this approach to circularly polarized signals, we replace 
F 2 a with F+ a + F 2 a in equation i@; this gives what we refer to as the circular 
PowerFlux method. The implementation of PowerFlux described in [H [3J uses the 
linear PowerFlux method with a search over discrete values of the polarization ip, 
together with the circular PowerFlux method, to search for elliptically polarized 
signals. 

3. Generalization to estimate A 2 + and A\ 

We now investigate generalizing the above derivation to estimate A 2 + and A 2 X 
simultaneously. The natural generalization of equation Q is to redefine g as 

_ ^ [P a - 0.5(A 2 + F 2 a + A 2 X F 2 Q )T SFT } 2 

a L si ■ (n 

Following the same minimization procedure, we obtain the equations 

dg = v [Pa - 0-HA 2 + F 2 a + A 2 x F 2 a )T SFT }F 2 a T SFT 

dg = v [Pg - Q-HA 2 + F 2 a + A 2 x F 2 a )T SFT }F 2 a T SFT 
dAl ^ si 

cx 

which we solve for and A 2 ^ to give 



+ V 

x v 



where 
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^ = E#E#- E^M • (12) 



We refer to this method as generalized PowerFlux I. A natural detection statistic 
would be A\ + A 2 , ; to evaluate 
the linear PowerFlux method. 1 
search over discrete values of ip. 



would be A\ + A 2 , ; to evaluate it, we must compute 5/2 as many summations as for 
the linear PowerFlux method. We must also, as for linear PowerFlux, still include a 
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4. Generalization to estimate A 2 ,, A\, and ip 

To further generalize the PowerFlux method to directly estimate ip, we note that the 
antenna patterns F+ and F x can be written in terms of two functions a and b, given 
in pQ, which are independent of ip: 

F + (tp,t) =sinC[a(i)cos2^ + 6(t)sin2^], (13) 
F x (ip,t) = sin cos 2<^ - a(t)sm2ip]. (14) 

where £ is the angle between the arms of the interferometer. The normalized signal 
power can therefore be written as 

= 0.5(A4 + Bb 2 a + Ca a b a )T SFT , (15) 

JSFT 

where the amplitudes A, B, and C are defined to be 

A = sin 2 C{A\ cos 2 2ijj + A\ sin 2 2^), (16) 
B = sin 2 C(A^ sin 2 2ip + A% cos 2 2ip), (17) 
C = sin 2 ({A 2 + - A\ )2 cos 2^ sin 2^. (18) 

Following the minimization procedure again, we redefine g to be 
^ [P a - 0.5{Aa 2 a + Bb 2 a + Ca a b a )T SFT } 2 

minimizing g with respect to A, B, and C gives 

dg [Pa - QHAal + Bb 2 a + Ca a b a )T SFT ]a 2 a T SFT 



OA ^ S 2 

a. " 

dg_ [P a - OMAaj + Bb 2 a + Ca a b a )T SFT ]b 2 a T SFT 

dB ^ S 2 



0, (20) 
0, (21) 



d£ _ _ sr [Pa - OMAaj + Bb 2 a + Ca a b a )T SFT }a a b a T SFT _ 

DC ~ ^ S 2 -U-(^) 

a 

Thus, the amplitudes A, B, and C can be found by inverting equations (|2"U)) - (l2"2"|) ; 
the amplitudes A 2 ^ and A\ and polarization angle ip are then found by inverting 
equations (fTI>)) - (fT8"|) . This method is referred to as generalized PowerFlux II. It involves 
computing 8/2 = 4 as many sums as the linear PowerFlux method; however, because 
we no longer need to search over discrete values of ip, the total computational cost of 
this method may overall be lower than the implementation of PowerFlux described in 
[3113]. An alternate method of estimating A 2 + , A\ , and ip, which uses the output of the 
linear PowerFlux method evaluated at several fixed values of "0, is given in Appendix 
A of [3J. 



5. Comparison of detection efficiencies 

We have determined the relative detection efficiencies (defined below) of the standard 
PowerFlux methods, and the generalizations presented above. The detection statistic 
used for each method is A 2 ^ + A\ , with A\ = in the case of the linear PowerFlux 
method. For simplicity, we did not restrict the estimations of A\ and A\ to the 
physical region A , > and A\ > 0, or check if the actual minimum value of g was on 
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Figure 1. Detection efficiency, determined to within 3% by simulated searches, 
versus the normalized injected amplitude, for 336 SFTs and signal declinations 
of zero (left) and 90° (right). The marked points are the results of the simulated 
searches and the solid lines are spline fits to these. 



the boundary of this region. An implementation of the StackSlide method described 
in [3], which uses a sum of the power as the detection statistic, was also included. 

We first performed searches on 10000 sets of SFTs containing only randomly 
generated Gaussian noise and obtained, for each method, a distribution of their 
detection statistic in the absence of a signal. From this distribution we determined the 
threshold for a 1% false alarm rate. We then performed searches on 3000 sets of SFTs, 
each containing noise plus a simulated signal from an isolated spinning source. In this 
case the strength of the signal can be given by the normalized injected amplitude 
M?sft/<S) 1/2 , where h is defined by PQ 

A + = 0.5/i (l + cos 2 i), (23) 
A x = h cost, (24) 

and where i is inclination angle between the spin axis of the source and the direction 
from the source to the detector, and S is the one-sided power spectral density of the 
Gaussian noise. For each of the 3000 simulated signals, cost was chosen at random 
from the range [—1,1], resulting in elliptically polarized signals; the polarization angle 
tp was also chosen at random from the range [— 7r/4, 7r/4]. 

Other parameters of the signal are its sky position (right ascension and 
declination) and its frequency at the Solar System Barycenter (SSB), taken as constant 
here. Since right ascension has minimal effect on the results this was set to zero and 
results were found for several values of the declination. For the detector configuration 
we used that of LIGO Hanford, as given in jT|. We also varied the number of SFTs. 

A real search requires a template bank of sky positions, frequencies (and frequency 
time derivatives, which we set to zero) at the SSB, and in some cases polarization 
angles. Here, for each simulated signal a single template was created with no mismatch 
in sky position. However, we included a random mismatch between the frequencies of 
the template and the signal at the SSB of up to half an SFT bin width and, for searches 
requiring a value for ip, a random mismatch in if) between the template and the signal 
of up to 7r/16. For each simulated source, the detection statistic was computed using 
the template sky position, SSB frequency and ip (if needed); thus included in the 
results are losses due to mismatch between the template and the signal. 
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Figure 2. Same as figured] for 336 (left) and 672 SFTs (right), and a signal 
declination of 45°. 

After searching the 3000 sets of SFTs we obtained, for each method, a distribution 
of their detection statistic in the presence of a signal of normalized amplitude 
^o(7sft/<S') 1 / 2 - The detection efficiency is then calculated as the fraction of the 
distribution of the detection statistic which falls above the 1% false alarm rate 
threshold. We then repeated the process of generating and searching 3000 sets of SFTs 
for increasing values of the normalized amplitude; this gives the detection efficiency 
as a function of h (T SFT / S) 1 / 2 . 

To verify our results, the authors independently wrote Matlab scripts to perform 
the above procedure. The scripts produced identical detection efficiency curves for 
the same input parameters, within the expected uncertainties based on the number of 
searches performed; this gives us confidence that our implementations are correct. 

Figures Q] and [2] show the detection efficiencies of the methods versus 
^o(7sft/<S') 1 ^ 2 , for a selection of values for the number of SFTs and declinations. As 
the source moves away from zero declination, the detection efficiency increases, and 
the differences between the methods becomes smaller. Also, the relative efficiencies of 
StackSlide, the generalized PowerFlux, and the linear PowerFlux methods can change 
with respect to each other, but note that circular PowerFlux remains typically the most 
efficient method. This is good news, since this method is already being used by the 
LSC. It also may not be surprising, since we note that the weights used in the Hough 
transform search in [3] were proportional to F? + F 2 , . These weights were shown to be 
optimal in an average sense in [3 [7], from which we see that A\F\ + A^F 2 . averaged 
over ijj can be factored as 0.5(F 2 + F 2 ){A 2 + + A 2 X ), i.e., as the antenna pattern for 
circular polarization times the sum of the squared amplitudes, even for an elliptically 
polarized signal. Since our simulations perform a Monte Carlo average over ip, this 
may explain why circular PowerFlux does so well. Finally, note that the detection 
efficiencies increase with the number of SFTs, as expected, but the relative efficiencies 
of the methods do not change significantly. 

Figure [3] shows, for a search of 3000 sets of SFTs injected with /i (T SFT /5) 1 / 2 = 4 
signals, the distribution of the injected versus the detected amplitudes estimated by 
the various methods. The mean value of the squared amplitudes estimated for noise 
alone was subtracted from the detected squared amplitudes; the injected and detected 
squared amplitudes are then normalized by Tsft/5*- The detected squared amplitudes 
are typically smaller than that injected due to loss of power due to mismatch between 
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Figure 3. Injected vs. detected (after removing the noise only mean) normalized 
squared amplitudes estimated by (from left to right, top to bottom) linear 
PowerFlux circular PowerFlux (A+ = A x ), generalized PowerFlux I (A+ 

and A x ), and generalized PowerFlux II (A + and A x ), for /io(Tg PT /S) 1 / 2 = 4. 
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the signal and the template used to find it, while the width of the distributions are 
due to noise. The generalized PowerFlux II method has marginally the narrowest 
such distribution. Figure 0] shows the injected versus the detected polarization angle 
estimated by generalized PowerFlux II; the estimation can be quite poor, for example 
when noise causes the detected ip to show up in the wrong quadrant. 



6. A PowerFlux maximum likelihood statistic, and future work 

In section [31 we gave equation for g, in analogy to the x 2 statistic. A direct 
comparison fails because P a is not a Gaussian distributed variable. However, if we 
expand the square in the numerator, and retain only the terms that depend on A + 
and A x , we can define a "PowerFlux maximum likelihood statistic", Q, by 

[(AlFf a + A\Fl a )T SFT P a - Q.2h{A\Fl a + A|igJ^y 
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where it is understood that A 2 ^ and A\ are chosen to minimize g. For example, we 
can substitute equation © for A\ and set A\ = 0, which gives 



Here, Q represents our definition of the "maximum likelihood statistic" for linear 
PowerFlux. It is similar to the standard linear PowerFlux statistic given in 
equation note, however, that the sum in the numerator is squared. Similar 
expressions for the circular and generalized PowerFlux methods, and indeed any 
method that computes either A 2 + or A\ , can be found using equation (|2"5j) . 

For future work we plan to investigate whether using the PowerFlux maximum 
likelihood statistic, Q, gives a better detection efficiency than the sum of the 
squared amplitudes, A 2 + + A\ , used in this paper. It would also be interesting to 
further understand why circular PowerFlux is so efficient (and under precisely which 
conditions this is so) , and study whether this method or another one is mathematically 
the optimal filter of SFT power in the Neyman-Pearson sense, i.e. that maximizes the 
detection efficiency for a fixed false alarm rate. 
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